global a0 a1
myfun3 = @(x) a0(3) + a0(2).*x + a0(1).*x.*x;
myfun4 = @(x) a1(3) + a1(2).*x + a1(1).*x.*x;
ta=[707.85 706.696 706.314 706.256 706.628 707.748 710.115];
me=[1051.51682334 977.396926963 935.398704788 927.021773958 948.615956238 1000.33563744 1080.0572534];

tb=ta-360; %shifting
tb=ref-568.39;
mw=me-568.39; %normalization
m=122:128;

a0=polyfit(m,tb',2);
a1=polyfit(m,mw,2);

figure(1);
hold on;
xlim([121.5 128.5])
xlabel('Mass (GeV)');
ylabel('-ln(L)');
scatter(m,tb,'vr','filled');
scatter(m,mw,'^','black','filled');
xpoints=121.75:0.1:128.25;
plot(xpoints,quadraticFunc(xpoints),'black');
plot(xpoints,quadraticFunc1(xpoints),'red');
legend('GPU','Original');

[min, minval] = fminunc(myfun4,125)





